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Multircsolution  analysis  of  SAR  Data 

Final  Report 

Robert  Hummel 
New  York  University 


1.  Overview 

The  “Multiresolution  Analysis  of  SAR  Data”  project  (Contract  F3361 5-89- 1087,  reference  (BAA) 
87-02-PMRE),  has  supported  research  work  in  five  areas: 

(1)  Geometric  hashing  foundations  and  theory  for  model-based  vision: 

(2)  Parallel  implementations  of  geometric  haslung; 

(3)  Motion  parameter  estimation  algorithms  for  egomotion; 

(4)  Uncertainty  reasoning  and  evidential  reasoning  studies;  and 

(5)  TTie  use  of  wavelets  for  SAR  image  construction. 

All  of  the  items  except  (3)  are  reasonably  part  of  the  original  SOW.  Item  (3)  developed  due  to  the  work 
of  V.  Sundareswaran,  who  was  supported  by  this  grant  after  I.  Rigoutsos  won  support  from  an  IBM  Fel¬ 
lowship. 

The  following  are  the  principal  advances  that  have  been  made  in  each  of  these  areas: 

(1)  Geometric  hashing  theory  can  now  be  viewed  as  a  Bayesian  approach  to  object  recognition.  False 
alarm  rates  can  be  greatly  reduced  by  using  certain  enhancements  and  modifications  developeu 
under  this  project. 

(2)  Geometric  hashing  algorithms  now  exist  for  the  Connection  Machine.  Reco^ition  of 
synthetically-pnaduced  dot  arrays  has  been  demonstrated  using  a  model  base  of  1024  objects.  The 
woik  represents  a  substantial  advance  over  existing  model-based  vision  capabilities. 

(3)  We  have  developed  new  algorithms  for  deteraiining  the  translation  and  rotation  of  a  sensor  given 
only  the  image  flow  field  data.  These  are  new  algorithms,  and  are  much  more  stable  than  existing 
computer  vision  algorithms  for  this  task.  The  work  has  been  accepted  for  publication  in  a  major 
journal.  The  algorithms  might  provide  independent  verification  of  gyroscopic  data,  or  might  be 
used  to  compute  relative  motion  with  respect  to  a  moving  scene  object,  or  may  be  useful  for 
motion-based  segmentation. 

(4)  Our  theories  explaining  the  Dempster/Shafer  calculus  and  developing  new  uncertainty  reasoning 
calculi  have  been  extended,  and  presented  at  a  conference,  and  were  incorporated  into  the  Bcye- 
sian  interpretation  of  geometric  hashing  (see  (1)  above). 

(5)  We  have  developed  a  “Wavelet  Slice  Theorem,”  in  several  different  versions,  any  of  which 
yields  an  alternate  approach  to  image  formation.  The  result  may  well  provide  a  more  stable 
aK>roach  to  image  formation  than  the  standard  Fourier-based  projection  slice  theorem,  since 
interpolation  of  unknown  spectral  values  is  better-founded.  Moreover,  we  have  made  consider¬ 
able  progress  on  an  alternative  approach  to  SAR  image  formation  that  avoids  the  computationally 
expensive  steps  of  rerasterizing  return  data,  and  Fourier  transforms,  using  instead  a  more  complex 
chirp  function.  The  method  suggests  a  simple  real-time  implementation  that  builds  up  the  final 
image  incrementally,  during  the  fly-by. 

The  list  of  significant  technical  reports  and  publications  produced  as  a  result  of  the  research  support  of 
this  grant  are  listed  below, 

[1]  R.  Hummel,  “Geometric  Hashing  Methods  for  Object  Recognition,”  NYU  Proposal  to  DARPA  in 
response  to  BAA  90-15,  submitted  October  1990,  reviewed  and  approved  (for  future  funding). 
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January,  1991. 

[2]  I.  Rigoutsos  and  R.  Hximmel,  “Several  Results  on  Affine  Geometric  Hashing,”  Proceedings  of  the 
8th  Israeli  Conference  on  Artificial  Intelligence  and  Computer  Vision,  pp.  1-12,  December  1991, 
Tel  Aviv,  Israel. 

[3]  I.  Rigoutsos  and  R.  Hummel,  “A  Bayesian  Approach  to  Model  Matching  with  Geometric  Hash¬ 
ing,’’  NYU  Technical  Report,  1991,  available  by  anon>’mous  fq?  from  cs.nyu.edu  under 
/pub/iocal/hummel/papers/bayeshashing. 

[4]  1.  Rigoutsos  and  R.  Hummel,  “Robust  Similarity  Invariant  Matching  in  the  Presence  of  Noise;  A 
Data  Parallel  Approach,’’  Proceedings  8th  Israeli  Conference  on  Artificial  Intelligence  and  Com¬ 
puter  Vision,  pp.  27-43,  December  1991,  Tel  Aviv,  Israel. 

[5]  Y.  Lamdan  and  H.  Wolfson,  “On  the  Error  Analysis  of  Geometric  Hashing,’’  IEEE  Computer 
Society  Conference  on  Computer  Vision  and  Pattern  Recognition,  CVPR91, 1991,  pp.  22-27. 

[6]  H.  Wolfson  and  Y.  Lamdan,  "Transformation  Invariant  Indexing,’’  in  Invariance  in  Computer 
Vision,  (Eds.  J.  Mundy  and  A.  Zisserman),  pp.  335-353,  MIT  Press,  1992. 

[7]  H.  Wolfson,  “Articulated  Object  Recognition,  or.  How  to  Generalize  the  Generalized  Hough 
Transform,”  (with  A.  Beinglass),  Proceeedings  of  the  1991  IEEE  Computer  Society  Conference  on 
Computer  Vision  and  Pattern  Recognition,  pp.  461-466,  June  1991. 

[8]  I.  Rigoutsos  and  R.  Hummel  “Massively  Parallel  Model  Matching:  Geometric  Hashing  on  the  Con¬ 
nection  Machine,”  IEEE  Computer,  pp.  33-42,  February  1992. 

(91  I.  Rigoutsos  and  R.  Hummel,  “Implementation  of  Geometric  Hashing  on  the  Connection 
Machine,”  Proceedings  IEEE  Workshop  on  Directions  in  Automated  CAD-based  Vision  Systems, 
June  1991,  Maui,  Hawaii, 

(101  1.  Rigoutsos  and  R.  Hummel,  “On  a  Parallel  Implementation  of  Geometric  Hashing  on  the  Connec¬ 
tion  Machine,”  CIMS  Technical  Report  #554,  April  1991  (Extended  version  of  the  previous  paper). 

[11]  R.  Hummel  and  V.  Sundareswaran,  “Motion  estimation  using  global  image  flow  velocities,” 
ARVO  Abstract,  1991. 

[12]  V.  Sundareswaran,  “Egomotion  from  global  flow  field  data,”  Proceedings  of  the  Workshop  on 
Image  Motion  Analysis,  Princeton,  N.J.,  October,  1991. 

[13]  R.  Hummel  and  V.  Sundareswaran,  “Motion  parameters  from  global  flow  field  data,”  To  Appear  in 
IEEE  Transactions  on  Pattern  Analysis  and  Machine  Intelligence,  May  1993. 

[14]  V.  Sundareswaran,  “A  fast  method  to  estimate  sensor  translation,”  in  Proceedings  of  the  Second 
European  Cortference  on  Computer  Vision,  Santa  Margherita,  Italy,  May  1992. 

[15]  R.  Hummel  and  L.  Manevitz,  “Statistical  approaches  to  fusion  with  uncertainty,”  Session  MB4, 
Proceedings  of  the  IEEE  International  Conference  on  Systems,  Man,  and  Cybernetics,  Charlotts- 
ville,  Virginia,  October,  1991. 

We  next  give  technical  details  on  the  progress  in  each  area. 

2.  Geometric  Hashing 

Geometric  hashing  is  an  approach  to  model-based  object  recognition  that  is  analogous  to  the  hash¬ 
ing  approach  to  data  retrieval  in  computer  science.  The  idea  is  that  collections  of  object  features  are 
encoded  in  such  a  way  that  recognition  of  a  few  features  may  enable  indexing  into  a  table  that  then  desig¬ 
nates  a  list  of  possible  objects.  The  intent  is  to  reduce  the  complexity  of  search  processing  by  allowing 
scene  features  to  index  into  the  model  database,  rather  than  positing  matches  between  scene  features  and 
model  features  and  then  searching  for  a  hypothesis  that  can  be  verified. 

The  initial  idea  is  due  to  Jack  Schwartz,  who  applied  the  idea  to  curve  matching  for  robotic  applica¬ 
tions.  In  this  form,  the  idea  was  related  to,  but  somewhat  dissimilar  from,  the  idea  of  R-tables  in  the 
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generalized  Hough  transform.  Subsequently.  Haim  Wolfson  and  his  students,  working  with  Professor 
Schwartz,  extended  the  idea  to  other  object-recognition  tasks,  including  3-D  curve  matching,  afhne- 
invariant  object  matching  based  on  point  sets,  affine-invariant  matching  based  on  curves  using  concavi¬ 
ties  as  primitive  features,  matching  using  convex  curve  portions  for  recognition,  and  matching  of  objects 
based  on  line  extracdons.  A  milestone  was  reached  with  the  publication  of  the  thesis  of  Y,  Lamdan,  at 
NYU,  who  studied  hashing  methods  for  3-D  object  recognition  from  2-D  projections.  His  thesis  is 
lecomraeikied  reading  for  an  understanding  of  the  basis  of  the  idea  and  the  many  variations  and  possibili¬ 
ties.  Tire  idea  has  been  taken  up  by  other  groups,  and  has  been  studied  by  researchers  at  the  University  of 
Washington,  MIT,  Syracuse  University,  and  in  Europe,  especially  at  INRIA-Rocqencourt. 

With  that  background,  and  support  from  AFAL,  we  have  subsequently  developed  the  theory  of 
geometric  hashing,  with  the  following  intents: 

(1)  We  have  analyzed  the  stability  and  false  alarm  rates,  both  experimentally  and  theoretically; 

(2)  Concentrating  on  rigid-transformation  and  similarity-transfonnation  invariance,  we  have  imple¬ 
mented  systems  using  thousands  of  models  of  dot  patterns; 

(3)  In  order  to  improve  the  efficiency  of  the  method,  we  have  developed  a  number  of  enhancements, 
including  (1)  Hash  table  equalization,  (2)  Improved  hash  functions  and  use  of  symmetries  in  the 
hash  table,  and  (3)  folding  of  the  hash  table. 

(4)  We  have  developed  a  Bayesian  interpretation  of  geometric  hashing  that  permits  a  well-founded 
basis  for  weighted  voting  in  the  evidence-gathering  phase  of  geometric  hashing. 

(5)  We  have  developed  the  application  of  geometric  hashing  to  articulated  object  recognition,  which 
should  facilitate  recognition  of  models  which  cen  vary  according  to  one  or  a  few  parameters. 

All  of  this  work  is  described  in  detail  in  the  publications  [1-7]  (from  the  list  of  references  in  Section 

1). 

The  application  of  geomeiric  hashing  to  ATR  problems,  both  with  SAR  data  and  IR  sensor  data,  is 
clear.  However,  we  have  not  pursued  these  issues  extensively  in  this  project,  because  the  development  of 
such  applications  is  critically  dependent  on  stable  feature  extraction  methods.  Such  algorithms  are  being 
extensively  developed  by  DoD  contractors,  but  require  large  databases  for  evaluation.  Much  of  our  work 
has  used  synthetic  data,  based  on  point  feature  sets,  with  noise  added,  in  order  to  evaluate  the  power  of 
the  indexing  methods.  Recently,  we  have  begun  to  experiment  with  recognition  of  military  aircraft  from 
silhouettes  and  extracted  edge  features.  However,  even  in  this  application,  our  data  comes  from  freely- 
available  pictures  and  sketches.  Rigoutsos  will  also  use  images  of  automobiles  in  his  continuing  thesis 
woik.  These  objects  are  somewhat  more  difficult  because  of  their  more  rotationally  symmetric  shapes. 
That  is,  the  military  aircraft  tend  to  be  elongated  with  distinguishing  marks.  To  date,  we  achieve 
extremely  good  results  for  similarity  recognition  of  military  aircraft  using  side-views,  employing  a  stan¬ 
dard  edge  detector  and  automatic  point  extraction. 

Most  of  the  developments  are  incremental,  with  the  possible  exception  of  the  Bayesian- 
interpretation  work,  which  might  be  considered  more  revolutionary.  Accordingly,  we  explain  the  Baye¬ 
sian  interpretation  in  general  terms  here. 

In  similarity-invariant  matching  of  point  sets  using  geometric  hashing  (which  we  use  here  as  an 
example),  pairs  of  points  are  chosen  from  the  scene  as  a  basis  pair.  For  each  basis  pair,  a  coordinate  sys¬ 
tem  is  defined,  and  the  other  scene  points  can  be  expressed  as  coordinate  pairs  in  any  given  coordinate 
system.  For  a  given  basis,  each  point  thereby  hashes  to  a  location  in  a  two-dimensional  hash  space.  In 
classical  geometric  hashing,  as  it  has  been  described  heretofore,  the  hash  location  is  regarded  as  landing 
within  a  bin,  and  each  entry  that  landed  within  that  bin  during  a  preprocessing  phase  receives  a  vote. 
However,  during  the  preprocessing  phase,  each  entry  is  associated  with  a  position  within  the  two- 
dimensional  hash  space.  Rather  than  regarding  the  entries  as  belonging  to  bins,  we  regard  them  as  fixed 
points  in  the  hash  space.  Accordingly,  when  a  point  hashes  into  the  space,  it  should  vote  for  nearby 
entries.  Further,  it  seems  reasonable  that  the  weighted  degree  of  vote  should  depend  on  the  distance 
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between  the  hash  point  and  the  entry.  The  larger  the  d-stance.  the  smaller  the  vole  should  be.  The  '■ate  of 
fall-off  of  the  vote  should  depend  on  the  expected  degree  of  noise  in  the  location  of  the  hash  points, 
which  in  turn  depends  on  the  expected  noise  in  the  placement  of  the  points  in  the  scene,  and  in  the  dis¬ 
tance  between  the  points  in  the  basis  pair.  When  one  applies  a  strict  Bayesian  analysis  to  the  .situation, 
the  result  is  a  precise  formula  for  the  weighted  contribution  that  a  hash  to  location  («,v)  should  make  to 
an  entry  for  model  k  using  basis  (i,y)  located  at  a  point  {x,y)  in  hash  space.  The  contribution  depends  on 
the  standard  deviation  of  the  noise  a  in  the  scene,  and  the  basis  separation  D.  In  this  case,  the  contribu¬ 
tion  is: 

Figure  1  shows  an  example  of  recognition  of  a  military  aircraft  from  a  database  of  1 5  aircraft.  The  scene 
comes  from  a  grayscale  image  of  one  of  the  aircraft  types,  but  is  not  taken  from  the  same  data  that  are 
used  to  construct  the  model  database.  The  points  in  the  scene  are  extracted  automatically  in  this 
example.  The  results  are  easily  extended  to  affine  invariant  matching. 

We  believe  that  this  work  forms  the  leading-edge  of  the  object  recognition  research  domain,  as  was 
evidenced  by  the  results  presented  at  die  1991  Workshop  on  Automated  CAD-based  vision  sy.stems.  The 
theoretical  advantages  of  the  geometric  hashing  method  are  clear,  and  the  issue  remains  as  to  verify  the 
applicability  through  experiments  and  through  refinements.  Most  other  researchers  are  working  with 
databases  consisting  of  one  or  two  objects,  and  most  work  remains  in  the  realm  of  geometric  reasoning 
and  hypothesis-and-verify.  Some  groups  have  begun  using  geometric  hashing  methods  and  their  variants, 
and  we  have  benefited  from  their  work,  but  the  NYU  work  remains  more  advanced  in  terms  of  both  appli¬ 
cations  and  generality. 

3.  Parallel  algorithms  for  geometric  hashing 

A  major  advantage  of  the  geometric  hashing  approach  to  object  recognition  is  its  parallelizability. 
However,  there  are  many  different  ways  to  parallelize  the  algorithm.  Chosing  an  appropriate  parallel 
implementation  can  be  very  subtle. 

Rigoutsos,  as  part  of  his  thesis  work,  has  implemented  a  number  of  versions  of  geometric  hashing 
for  execution  on  the  Thinking  Machines  Inc.  Connection  Machine™.  We  have  gained  access  to  various 
Connection  Machines  by  making  use  of  the  DARPA  “Network  Server”  program.  He  is  currently  com¬ 
pleting  an  updated  implementation  of  one  of  the  versions  which  will  support  the  weighted  voting  as 
described  in  the  Bayesiart  imeipieuuou  above,  and  will  be  used  fvi  th '  rapid  recogniiio”  of  military  air¬ 
craft  and  automobiles. 

As  is  very  typical,  there  is  an  interesting  and  subtle  parallel  approach  which  has  good  asymptotics, 
and  there  is  a  much  simpler  and  less  sophisticated  algorithm  that  is  likely  to  present  the  more  practical 
approach  in  general.  The  sophisticated  approach  uses  a  “connectionist”  view  of  geometric  hashing, 
whereas  the  simpler  version  is  based  on  a  'hash  location  broadcast”  model.  The  two  approaches  are  best 
understood  from  the  paper  that  will  appear  in  IEEE  Computer  (Reference  [8]).  This  paper,  which  will  be 
a  very  high-profile  presentation  of  the  geometric  hashing  paradigm  and  the  parallelizability  of  the 
methods,  is  a  relatively  short  and  “popular”  description,  intended  for  a  very  broad  audience  (all  members 
of  the  IEEE  Computer  Society  receive  the  journal).  References  [9]  and  [10]  further  describe  aspects  of 
the  parallel  implementation. 

The  fundamental  choice  in  our  parallelizations  of  geometric  hashing  is  to  continue  to  serialize  over 
choices  of  basis  pairs  in  the  scene,  but  to  maintain  the  entire  hash  table  data  base  in  the  memory  of  the 
parallel  machine,  and  to  make  sure  that  the  search  over  model/basis  pairs  is  done  in  parallel.  This  choice 
is  motivated  by  the  fact  that  the  advantage  of  geometric  hashing  is  that  it  makes  the  search  over 
modelAtasis  pairs  fast,  whereas  it  retains  the  need  to  search  over  the  image  for  a  basis  set.  However  logi¬ 
cal  this  choice  seems,  all  other  parallel  implementations  of  geometric  hashing  of  which  we  are  aware 
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Figure  1.  Recognition  of  a  military  aircraft.  The  database  consists  of  15  military  aircraft  types,  where 
points  were  selected  manually  from  scanned  images  of  side-view  drawings  of  each  aircraft.  The  scene  is 
a  photograph  of  one  of  the  aircraft  types,  in  flight.  The  photograph  comes  from  a  different  book  than  the 
sketches  that  are  used  for  the  model  base.  Edges  and  points  are  extracted  automatically  from  the  scene, 
and  the  Bayesian  version  of  the  geometric  hashing  method  is  used  for  recognition.  The  nine  top 
model/basis  pairs  are  shown;  the  one  receiving  the  greatest  weight  is  the  correct  recognition  of  the  air¬ 
craft.  Bars  over  the  sketches  indicate  the  degree  of  support  for  each  model/basis. 


(including  the  Martin  Marietta  Denver  approach  as  of  last  summer),  parallelize  over  a  different  data  set, 

Allow  us  to  explain  in  somewhat  greater  detail  the  broadcast  approach.  We  have,  from  the  prepro¬ 
cessing  phase,  a  large  number  of  records  of  the  form  (Mic,B,l,x,y,C),  where  Mt  is  a  model  number,  S  is  a 
basis,  within  the  model,  /  is  another  point  in  the  model  (not  in  the  basis),  {x,y)  is  the  location  to  which  the 
point  hashes  under  the  basis  B,  and  C  is  a  covariance  matrix  which,  based  on  expected  noise  of  placement 
of  points  in  the  scene,  predicts  the  covariance  of  the  distribution  of  the  location  of  the  hash  value  of  the 
point  /  about  tte  location  (x.y)  under  the  basis  B.  The  hash  table  is  no  longer  organized  as  a  table,  but 
rather  as  simply  a  large  collection  of  such  records. 
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During  the  recognition  phase,  a  basis  is  chosen  and  fixed.  We  operate  sequentially  over  the  other 
points  in  the  scene.  It  is  this  degree  of  scrialism  that  one  would  like  'o  parallelize,  and  leads  to  thi  more 
complicated  connectionism  algorithm  ;»i  the  hroadcast  algorithm,  each  point  in  the  scene  successiycly  is 
hashed  to  a  location  (u.v)  in  the  ^  '■'t  tabic,  and  the  location  (m.v)  is  broadcast  to  all  the  records  stored  m 
the  hash  table.  (Some  saving-  at  be  obtained  by  broadcasting  to  only  those  records  that  lie  near  (u,v ), 
which  means  that  the  records  must  be  organized  in  a  suitable  manner  in  advance.)  Each  record  computes 
the  distance  between  its  (x.y)  and  the  location  (a.v),  and  uses  the  predicted  covariance  C  to  compute  the 
degree  of  contrib'..uon  that  the  scene  point  should  lend  to  the  as.sociated  information  Votes  arc 

accumulated  tnd  after  cycling  through  all  appropriate  scene  points,  votes  are  summed  for  modci/basi.s 
pairs  (A/i  ‘  These  form  the  desired  evidence  values. 

Figure  2  shows  an  example  of  a  recognition  that  takes  place  using  the  broadcast  algorithm,  as 
implemented  on  a  1 6K-processor  Connection  Machine. 

Many  variations  are  possible.  However,  as  stated,  the  algorithm  is  able  lo  work  efficiently  and 
rapidly  on  an  8K-processor  CM-2  with  16  models,  with  approximately  20  points  in  each  model,  using 
scenes  of  roughly  100  points.  Currently,  each  probe  takes  about  2  seconds,  and  typically  .10  probes  arc 
needed  (i.e,,  30  choices  of  basis  pairs)  in  order  to  find  the  model.  Approximately  linear  speedup  can  be 
achieved  on  a  larger  Connection  Machine,  (e.g.,  a  64K  model),  either  a  much  more  complicated  algo¬ 
rithm  must  be  used,  or  the  problem  must  be  scaled  up.  The  current  implementation  on  the  HK-procc.ssor 
machine  can  handle  32  models  in  the  database  without  any  degradation  in  processing  time. 
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Figure  2.  An  example  recognition  using  the  broadcast  algonihm  on  a  Connection  Machine.  There  arc 
1024  models  in  the  database.  Each  model  consists  of  16  points.  The  scene  contains  200  points.  Approxi¬ 
mate  0.8  second  was  required  for  each  probe  on  a  I6K-processor  CM.  By  using  symmeuics  in  the  hash 
table,  this  time  could  easily  be  halved.  A  roughly  fourfold  speedup  would  be  expected  by  using  a  64K- 
proccssor  CM.  For  this  example,  150  probes  were  required  before  the  model  was  found.  However,  fewer 
probes  arc  typically  required  when  the  data  is  not  based  on  synthetic  models. 


Using  the  connectionist  parallel  approach,  and  databases  of  512  objects,  with  16  poiiiis  in  each 
object,  and  scenes  of  200  points,  we  have  achieved  probe  times  of  under  200  milliseconds  on  a  22 
processor  CM-2. 

In  any  case,  we  are  able  to  perform  recognition  of  models  in  less  than  a  minute,  wiiii  dozens  ol 
models  in  the  database,  which  is  considerably  beyond  the  state  of  the  art  in  model-based  object  recogni¬ 
tion,  and  we  are  satisfied  that  real-time  recognition  using  these  methods  is  possible  with  the  appropriate 
parallel  hardware. 

4.  Motion  parameter  estimation  from  optical  flow 

This  work  is  considerably  unrelated  to  the  Statement  of  Work  of  the  contract,  and  was  only  lightly 
supported  by  the  contract.  However,  because  student  V.  Sundareswaran  became  interested  in  analyzing 
sensor  motion,  and  because  we  wished  to  develop  global  methods  that  support  a  greater  degree  of  stabil¬ 
ity  than  existing  algorithms,  we  pursued  the  topic  of  motion  parameter  estimation. 

This  is  a  problem  that  has  concerned  the  lU  community  for  more  than  a  decade.  Immense  amounts 
of  effort  have  been  expended  on  motion  analysis.  Ttere  are  two  aspects  the  problem.  First,  giver,  an 
image  sequence,  we  typically  wish  to  understand  the  image  motion  from  p  xcl  to  pixel.  The  result  is  a 
vector  field,  called  tlie  image  flow  field  or  the  optical  fiow  field.  The  idea  ij  that  every  pixel  points  to  the 
corresponding  pixel  in  the  next  frame,  which  moves  according  to  the  motion  of  the  sensor  (and  pciien- 
tially  the  motion  of  the  object  imaged  at  the  pixel).  Second,  there  is  the  problem  of  determining  the 
motion  parameters  from  the  optical  flow  field.  This  is  an  inverse  problem,  in  that  the  optical  flow  is 
determined  by  the  instantaneous  translation  and  rotation  of  the  sensor  (and  objects).  In  addition,  there  is  a 
dependence  on  the  distance  from  the  object  to  the  sensor.  Accordingly,  one  can  hope  to  reconstruct  the 
motion  (translation  and  rotation)  parameters,  and  the  depths  to  the  objects. 

What  has  interested  most  researchers  in  this  topic  is  the  depth  dependence.  1  heoreiically.  one  can 
determine  the  optical  flow,  analyze  the  result  mathematically,  and  determine  the  parameters  of  motion  of 
the  sensor,  the  depth  to  scene  points,  and  identify  moving  objects.  The  last  two  objectives  arc  useful  for 
image  segmentation  and  for  object  identification,  both  of  which  are  fundamental  to  ATR  applications. 

Alas,  all  of  the  algorithms  are  terribly  unstable.  There  is  some  hope  of  being  able  to  segment  out 
moving  objects,  using  work  recently  developed  at  Samoff  Research  Labs,  but  the  task  of  extracting 
depths  from  image  motion  flow  alone,  using  any  of  the  existing  algoriihn  s,  is  completely  hopeless.  The 
lack  of  stability  is  not  usually  admitted  or  discussed  prominently,  but  is  evident  from  the  methods.  The 
difficulty  is  that  the  depths  must  be  determined  from  local  data,  whereas  the  optical  flow  field  can  only  be 
determined  approximately,  sparsely,  (for  example  where  there  are  moving  point  features'),  and  only  'Aiih 
noise.  If  one  is  to  have  any  hope  of  making  use  of  optical  flow  field  data,  then  the  algorithm  must  be  glo 
bal. 

Globally,  only  the  motion  parameters  can  be  determined.  However,  once  motion  parameters  for  the 
sensor  are  detennined,  it  should  be  considerably  easier  to  extract  moving  objects  (providing  they  arc 
small,  and  thus  have  not  significantly  degraded  the  motion  parameter  estimation),  and  the  depths 
Accordingly,  we  have  developed  global  methods  for  motion  parameter  estimation  from  optical  flow- 
fields. 

We  developed  two  methods  for  the  determination  of  the  parameters  of  motion  of  a  sensor  given  the 
vector  flow  field  induced  by  an  imaging  system  governed  by  a  perspective  transformation  of  a  ngid  scene. 
We  assume  that  the  flow  field  V  =  (u(x,y),v(x,y))  is  given.  Both  algorithms  are  new,  and  both  arc 
c’^lr^ely  simple.  The  first  algorithm  is  called  the  flow  circulation  algorithm,  and  determines  the  rota¬ 
tion^  parameteis.  It  uses  the  curl  of  the  flow  field,  curKV).  We  do  not  present  details  here,  but  wc  show 
in  our  papers  that  under  many  conditions,  the  curl  is  approximately  a  linear  function  of  the  form  g(x.,y)  - 
ax  +  by  +  c  (where  x  and  y  are  image  coordinate  values).  The  coefficients  of  the  linear  function,  a.  h.  and 
c,  which  may  be  determined  by  simple  regression,  are  proportional  to  rhe  desired  rotational  parameters  of 
motion.  Circulation  values  may  be  used  in  place  of  curl  v^ues,  resulting  in  less  noise. 
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The  second  algorithm,  the  FOE  search  algoriihm,  determines  the  translaiiona]  para  -icrs  of  the 
motion  independently  of  the  first  algorithm.  Tliis  algorithm  extends  a  recently  introduced  mcihtxl  of 
Hecger  and  Jep.son,  giving  a  method  for  searching  for  the  image  focus  of  expansion  SpecilicaHy,  for 
every  location  (Xo.yo)  in  the  image  plane,  we  compute  the  function  u  (->+>o)  +  v{x~xq),  which  i.s  the 
component  of  flow  orthogonal  to  the  radial  direction  scaled  by  the  distance  from  the  point.  When  cto  To) 
is  located  at  the  focus  of  expansion,  this  function  will  be  a  quadratic  polynomial  (of  a  special  form*).  The 
algorithm  determines  for  every  (j:o->'o)  whether  the  computed  function  is  a  quadratic  polynomial,  in  order 
to  locate  the  focus  of  expansion.  We  suggest  several  candidate  methods  for  determining  when  the  func¬ 
tion  has  the  appropriate  form:  one  method  involves  filtering  the  function  by  a  collection  of  of  circular- 
surround  zero-mean  receptive  fields.  The  other  methods  project  the  function  onto  a  linear  space  of  qua¬ 
dratic  polynomials,  and  measures  the  distance  between  the  two  functions.  Each  suggested  method  can  be 
formulated  as  the  evaluation  of  a  quadratic  functional  at  each  candidate  point,  and  the  first  two  methods 
permit  a  quadratic  search  for  the  focus  of  expansion,  yielding  an  especially  rapid  search. 

It  is  this  final  point,  the  quadratic  nature  of  the  error  surface,  that  is  the  most  startling  result  of  this 
work  The  fact  that  there  exists  a  global  method  for  the  determination  of  motion  parameter  is  not  a 
surprise  —  it  is  only  a  surprise  that  these  methods  had  not  been  previously  developed  and  exploited. 
How  ever,  the  simplicity  and  reduction  in  cost  due  to  projecting  onto  a  space  of  quadratic  polynomials  are 
great  advantages  of  the  method,  and  argue  strongly  for  its  viability  as  part  of  a  motion  analysis  system. 

.\s  an  example,  we  show  a  contour  plot  (together  with  gradient  vectors)  for  the  polynomial  projec¬ 
tion  method  for  the  FOE  .search  algoriihm  applied  to  a  synthetic  flow  field  (.see  Figure  3), 

Either  algorithm  may  be  used  in  isolation  to  obtain  some  parameters,  from  which  the  c.bcr  parame¬ 
ters  may  then  be  easily  determined.  However,  the  two  algorithms  in  combination  provide  multiple  path¬ 
ways  for  independent  dcteimination  of  motion  parameters.  In  the  paper  accepted  for  publication  in  PAMI 
(reference  [  1 3]).  we  demonstrate  the  viability  of  both  approaches  with  a  number  of  calculations  for  lypi- 
caJ  situations,  and  a  number  of  experiments  using  synthetically-generated  motion  data.  Presentation  of 
other  details  of  the  methods  are  given  in  references  (11,12,14].  Sundareswaran  has  since  applied  the 
methods  to  real  motion  data  sequences,  and  has  investigated  sources  of  error  and  noise  in  the  system.  We 
are  currently  in  touch  with  NASA-Ames  researchers  who  will  supply  us  with  additional  data. 

Follow-up  work  would  anempi  to  make  use  of  the  motion  parameter  estimation  methods  to  locale 
moving  objects,  or  to  estimate  range  values,  so  as  to  determine  (for  example)  ground  plane  orientation, 

5.  Uncertainly  reasoning 

Continuing  a  longstanding  interest  in  uncertainty  reasoning,  Professor  Hummel  has  devoted  some  of 
the  ume  supponed  by  the  AFAL  contract  to  pursue  theories  of  information  fusion  and  uncertainty  calculi. 

Although  only  one  conference  paper  was  presented  (reference  [15]),  a  major  journal  article  is  in 
preparation,  and  will  be  submiiicd  soon.  All  the  work  is  joint  with  PTofcs.sor  Larry  Manevilz.,  at  Haifa 
University. 

Rca.soning  with  uncertainty  is  a  field  with  many  different  approaches  and  viewpoints,  with  impor¬ 
tant  applicatioas  to  .scn.sor  design  and  autonomous  system  development.  We  attempt  to  unify  some  of  the 
different  approaches  by  introducing  a  common  philosophical  framework  under  which  different  calculi 
may  N:  developed.  Each  calculus  reflects  different  design  choices  compatible  with  the  philosophical 
tenets.  The  tenets  postulate  that  uncertainty  in  representations  can  be  viewed  as  a  degree  of  dispersion  of 
opinions,  and  that  the  space  of  opinions  operates  as  a  separate  sample  space  distinct  from  the  underlying 
sample  space  on  which  probabilities  arc  normally  defined.  Different  calculi  result  when  choices  are  made 
for  the  representation  of  the  opinions,  the  method  for  combining  opinions,  the  method  for  juxtaposing 
multiple  sets  of  opinions,  and  the  way  of  measuring  the  spread  in  the  opinions. 

Many  sensor  fusion  systems  make  use  of  an  evidential  reasoning  system,  where  evidence  is  com¬ 
bined  With  current  measuremenLs  m  order  to  maintain  states  of  belief  and  confidence  in  a  set  of 
hypotheses.  These  systems  arc  all  motivated  by  the  fact  that  the  sensor  is  supposed  to  provide  more  than 


Figure  3.  Contour  and  gradient  plot  of  the  error  function  £(j;o,>o)  obtained  using  the  quadratic  polyno¬ 
mial  i»ojection  method  fa*  the  image  flow  field  induced  by  a  translational  velocity  of  (5,2,20)  and  rota¬ 
tional  velocity  of  (0.1,0.025,-0.05)  when  imaging  a  corridor  (four  side  walls  and  a  back  wall).  The  focus 
of  expansion  is  correctly  located  as  the  zero  of  the  error  surface.  It  can  be  seen  that  the  error  surface  is 
quadratic,  and  thus  the  zero  can  be  found  easily. 


simple  measurements,  but  also  confidence  levels  in  the  measurements,  and  ultimately  confidence  levels  in 
propositions  and  hypotheses  developed  from  those  measurements.  The  fundamental  concepts  always 
involve  quantities  related  to  the  degree  of  validity  of  a  proposition,  such  as  a  probability,  and  other  quan¬ 
tities  related  to  the  degree  of  certainty  in  the  assertion  of  the  degree  of  belief.  Various  calculi  arc  used  for 
representing  these  concepts  and  performing  the  calculations,  including  Bayesian  networks,  fuzzy  logic. 
Kalman  filtering,  and  the  Dempster/Shafer  theory  of  evidence.  Each  calculus  has  certain  theoretical 
underpinnings,  although  a  universally  accepted  methodology  is  still  lacking. 

Our  approach,  and  the  unified  treatment  of  calculi,  is  founded  on  certain  fixed  tenets.  The  tenets 
are: 

(1)  That  uncertainty  can  be  represented  by  a  distribution  of  opinions,  whereas  certainty  is  represented 
by  unanimity  of  multiple  opinions. 

(2)  By  an  “opinion,”  we  mean  an  estimate  of  a  quantity  that  is  functionally  related  to  a  (frequency- 
I  based)  probability,  or  is  a  subjective  estimate  of  a  conditional  probability,  or  is  an  estimate  of  a 

well-defined  quantity  representing  the  likelihood  of  a  given  proposition  based  upon  given  evi¬ 
dence. 

(3)  When  opinions  arc  combined  in  order  to  make  estimates  that  are  conditioned  on  combinations  of 
evidence,  a  precise  and  well-defined  combination  formula  should  be  used.  In  the  case  that  the 
opinions  are  estimates  of  probabilities,  Bayes’  rule  should  be  used  to  combine  pairs  of  opinions  to 
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yield  a  new  opinion. 

The  propositions  can  relate  to  measurement  values  obtained  from  different  sensors,  or  can  be  propositions 
that  are  developed  and  based  upon  those  measurements  (such  as  the  presence  or  absence  of  a  target  class). 
For  ATR  applications,  we  might  be  combining  information  from  different  modalities,  or  from  different 
target  recognition  algorithms,  or  from  different  environment  variables  together  with  feature  data. 

A  variety  of  different  calculi  are  obtainable  from  these  philosophical  foundations,  depending  upon 
the  values  that  the  opinions  are  supposed  to  represent,  and  depending  upon  the  assumptions  used  in  the 
updating  process.  Each  design  choice  leads  to  a  different  calculus,  giving  a  Chinese  menu  of  uncertainty 
calculi.  Different  applications  will  require  different  design  choices.  One  of  the  calculi  that  arises  is  pre¬ 
cisely  the  Dempster/Shafer  calculus.  Another  is  the  simplified  Kalman  filter. 


We  present  here  the  formulas  for  yet  another  calculus  that  arises  in  our  uniform  framework.  We 
envision  a  collection  of  “experts”  £,  wi^  each  expert  w  in  £  expressing  an  opinion  Xa,(X)  for  every  label 
Xe  A.  Suppose  that  the  opinions  are  represented  by  the  logarithms  of  values  determined  by  the  probabili¬ 
ties  for  a  particular  proposition.  The  values,  however,  are  not  simply  the  logarithms  of  estimated  condi¬ 
tional  probabilities  based  on  the  known  observations,  but  rather  the  log  of  the  ratio  of  the  conditional  pro¬ 
bability  and  the  prior  probability.  Specifically,  we  set: 


x^iX)  =  log; 


Prob(X|  Information) 
Prob(X) 


where  the  “Information”  is  the  information  shared  by  the  experts  in  £,  and  the  denominator  is  the  prior 
probability  of  label  X  in  the  absence  of  any  information.  Note  that  these  probabilities  are  defined  over  the 
usual  sample  space  of  problem  instances,  and  not  the  set  of  experts.  The  equations  are  approximate, 
because  each  expert  is  making  an  estimate.  These  values  are  the  representation  suggested  by  Chamiak  for 
probabilistic  reasoning.  Each  expert  will  have  a  different  estimate  of  this  log-ratio,  and  the  statistics  that 
we  propose  to  maintain  are  the  mean  and  variance  of  these  values.  Thus  if  an  expert  regards  a  proposition 
as  being  4  times  as  likely  due  to  the  given  measurements  as  opposed  to  its  probability  in  the  absence  of 
information,  then  the  expert  contributes  log(4)  as  the  opinion  in  the  set  of  opinions,  from  which  we  meas¬ 
ure  the  mean  and  variance.  Note  that  if  the  information  has  no  influence  on  the  prior  probability,  accord¬ 
ing  to  an  expert,  then  the  expert’s  opinion  will  be  log(l).  i.e.,  zero. 

The  Slate  of  the  system  is  thus  represented  by  two  vectors: 


^(X)  =  Avg(x„(X)). 
(0€£ 


and: 


a(X)  = 


AvB(xo>(X)-p,(X))^ 


1/2 


Unlike  the  Dempster/Shafer  representation,  which  requires  2^  values  for  the  specification  of  the  slate,  this 
formulation  requires  only  2N  values,  where  N  is  the  number  of  possible  labels  (i.e.,  N  =  #A). 


Now,  if  two  such  sets  of  opinions  are  to  be  combined,  we  chose  to  take  the  set  product  of  the  opin¬ 
ions.  If  a  composite  opinion  is  to  be  formed  from  two  individual  opinions,  and  if  we  may  suppose  a  con¬ 
ditional  independence  between  (he  information  sources  on  which  the  two  experts  are  basing  their  opin¬ 
ions,  then  it  can  be  shown  that  modulo  a  uniform  additive  constant,  the  two  opinions  may  be  summed. 
This  comes  from  Bayes  formula,  using  conditional  independence  of  the  information  sources,  yielding: 

Prob(X|  InfO|,Info2)  Prob(X|  Infoi)  Prob(X|  Info2) 

Prob(X)  ^  Prob(X)  Prob(X) 

The  two  sides  are  equal,  except  for  a  proportionality  constant,  which  is  independent  of  X.  Taking  the  log 
of  both  sides,  wc  sec  that  it  is  logical  to  set: 
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a)+x^{k). 

The  proportionality  constant  has  been  dropped,  which  means  that  the  opinions  can  be  off  by  a  constant 
additive  amount  (but  the  same  constant  for  all  the  labels  X).  This  skewing,  it  turns  out,  is  unimportant, 
siiKe  we  are  only  concerned  in  the  relative  size  of  the  components  over  different  labels  X.  This  then  is 
the  updating  method,  and  we  can  see  that  independence  (conditioned  on  every  label)  is  required  for  the 
informatioiL  A  mathematical  statement  of  the  independence  assumptions  says  that 

Prob(Infoi  I  Info2,  X)  =  Piob(lnfoi  I  X) 

for  all  X  in  A.  These  are  strong  requirements,  but  potentially  valid  in  some  circumstances. 

Using  the  product  formulation  for  obtaining  the  the  combination  set  of  opinions,  wc  find  that  the 
following  formulas  hold.  The  mean  and  standard  deviation  of  the  opinions  of  the  experts  in  set  E  j  are 
denoted  by  (pi  .Oi),  and  the  opinions  of  the  experts  in  set  E2  give  rise  to  the  state  (P2.O2)  Then 

^(X)  =  ^,(X)  +  ^t2(X), 

a(X)  =  j^Oi(X)  +  02(X) 

That  is,  the  resulting  calculus  can  be  specified  by  the  statement  that  mean  log-probability  opinions  should 
be  added,  and  variances  also  add.  We  look  for  a  situation  where  the  resulting  mean  opinion  is  cither  ver\ 
large  positive,  or  very  large  negative,  with  a  relatively  tight  variance  (much  less  than  the  magnitude  of  the 
opinion)  in  order  to  conclude  that  the  corresponding  proposition  is  tme  or  false. 

6.  Wavelet  slice  theorem  and  variants 

Our  work  on  wavelet  methods  for  SAR  image  fonnation  has,  to  date,  been  mostly  based  on 
mathematical  explorations.  The  work  was  begun  by  postdoctoral  student  Jacques  Froment,  and  has  been 
taken  up  by  mathematics  graduate  student  Jordan  Mann.  Wc  announce  here  principal  results.  The  first 
two  results  are  termed  “Wavelet  Slice  Theorems,”  and  they  form  wavelet  analogs  of  the  "Projection 
Slice  Theorem,”  which  forms  the  basis  of  the  SAR  technology  (as  well  as  tomography  technology).  The 
last  two  results  are  the  most  recent  work,  and  suggest  exciting  new  methods  for  practical  SAR  image  for¬ 
mation  that  might  provide  practical  advantages  over  current  methods. 

The  problem  of  reconstructing  a  function  of  n  variables  from  its  integrals  along  hypcrplanes  was 
first  studied  by  Radon.  The  problem,  its  two-dimensional  version  in  particular,  has  found  applications  in 
such  diverse  fields  as  radio  astronomy,  molecular  biology,  and  computer- assisted  tomography.  More 
recently  the  foundational  mathematics  has  been  applied  to  synthetic  aperture  radar  (SAR). 

We  have  proposed  several  new  formulae  which  hold  particular  promise  for  SAR.  since  they  take 
advantage  of  the  convolution  inherent  in  the  transmission  and  reception  of  the  radar  signal.  The  first  two 
formulae  are  based  on  arbitrary  wavelets  and  thus  open  the  possibility  of  determining,  analytically  or 
empirically,  the  wavelet  or  wavelets  most  effective  for  a  given  application.  Notably,  finely  tuned  orthogo¬ 
nal  wavelets  are  not  needed  in  our  methods. 

Since  these  ideas  are  as  yet  unpublished,  but  may  lead  to  exciting  new  methods  for  SAR  image  for¬ 
mation,  we  present  a  mathematical  summary  of  our  work  below. 

We  consider  a  function  /(x,y),  and  define  Fe/(r)  to  be  the  Radon  transform  of  the  function  at  angle 
0,  given  by 

f’e/  (t)=|2/(icos0-^sin0,rsin0+Tcos0)dx. 

For  .my  particular  t,  Pef(t)  is  the  integral  of  the  function  along  a  specified  line  in  the  coordinate  plane. 

Given  a  function  fix)  of  one  variable  and  a  1-D  wavelet  \|r(x),  we  define  WJiu),  the  wavelet 
tran.sform  of /at  location  u  and  scale  s  by: 


the  convolution  of  the  signal  / and  the  wavelet,  and  where; 

=  and 


Given  a  function  / ix,y)  of  two  variables  and  a  2-D  wavelet  ¥U<y)>  wavelet  transform  W^f  (u.v)  of  / 
at  location  (u.v)  and  at  the  scale  s,  is  defined  by: 

W'i/(w.vK/‘*¥^)(«.v). 

where  ¥j(^->)  and  ¥  defined  similarly  as  above,  but  for  two  variables. 

In  spotlight  mode  SAR,  if  we  “chirp”  a  function  h(t)  in  a  direction  6  toward  an  image  g(x,y),  the 
return  signal,  in  a  certain  sense,  is  essentially  h*P^g(t).  Accordingly,  in  the  foDowing  results,  we  con¬ 
sider  to  be  the  function  for  reconstruction. 

Result  1.  We  have  the  following  formula: 

c  f  ([/(e))-l^dsde = (m(ur)*g)(x,y) 

*  °  '^s  I— — - 

where  C  is  a  multiplicative  constant,  m(i/r)(x,y)  is  the  function  whose  Fourier  transform  is  l/'^^^+ri^ . 
and  6f(0)=jrcos0-fysin0. 

Since  the  convolution  on  the  right  side  of  the  equation  can  be  inverted,  we  can  recover  g.  In  fact,  we 
have; 

g  (x,y)  =  cQ^H'W.PQg  (V  (0))^  ds  d6. 

where  H'  denotes  differentiation  followed  by  the  Hilbert  trarisform. 

The  significance  of  this  result  is  that  we  may  use  any  1-D  wavelet  to  construct  W^P^g  by  chirping  a 
collection  of  different  signals  at  varying  scales  at  each  angle.  The  function  g  may  then  be  recovered 
through  the  indicated  formula,  which  involves  a  Hilbert  transforrri. 

Compared  with  typical  SAR  analysis,  a  rasterization  onto  a  coordinate  grid  is  not  necessary,  since 
all  the  processing  is  1-D  before  integrating  over  space.  However,  like  regular  SAR  processing,  FFTs  (or 
this  case,  Hilbert  transforms)  are  necessary. 

Result  2.  Let  ¥  be  any  2-D  wavelet.  By  the  wavelet  transform  inversion  formula: 

where  Cy  is  a  constant  associated  with  y.  writing  the  convolution  on  the  right  side  as  the  inverse 
Fourier  Transform  of  its  own  Fourier  transform,  and  the  convolution  theorem,  we  find  the  the  right  side  of 
(1)  is  equal  to; 

The  integral  over  d^d%  may  be  easily  transformed  into  a  polar  coordinate  integral  Ip/dpdB. 

By  making  the  chmge  of  variables,  and  regrouping  terms,  (1)  may  be  rewritten  as  follows.  Let; 

Cs  e(«)=-7==^  [  ¥j(pcos0.psin0)  \|r,(pcosO,psin0)/p/e‘P“  dp. 

V27t  ^ 

Then  (1)  becomes,  after  some  algebra: 

Wc  see  that  g  can  be  rcconstaictcd  using  an  extremely  simple  formula,  providing  we  have  available  the 
return  signals  from  a  spotlight  mode  radar  that  has  emitted  a  collection  of  Cj «  signals.  This  leads  us  to 
investigate  properties  of  C,  e 
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Suppose  we  set  Ceiu)  to  be  equal  to  C,.e(“)  when  s  =  1 .  Then  is  is  not  hard  to  show  that: 

C,_q(u)  =  C oisu).  and 
C  e(-u)  =  C  q(u). 

In  the  case  when  y  is  real  (and  real-valued  wavelets  are  the  normal  situation),  we  have; 

C'e(«)  =\~f  /V(pcose.psin0)/^ pcos(pu) Jp. 

’  It  ^0 

Basically,  C,_e(M)  is  a  function  whose  Fourier  Transfom  is  p  /\|r,/^(pcos6,psin0)  as  a  function  of  p. 
Accordingly,  we  expect  that  for  a  suitable  choice  of  y,  the  function  Q  q  wiU  be  smooth  and  essentially 
compactly  supported,  and  thus  a  good  choice  for  a  chirp  signal. 

The  problem  with  Result  2  is  that  an  entire  collection  of  different  scales  are  needed  at  each  angle. 
The  same  situation  occurs,  of  course,  with  standard  SAR.  The  issue  is  whether  there  exists  a  method  to 
chirp  a  range  of  different  frequencies,  and  sort  out  the  necessary  coefficients,  as  is  done  with  standard 
SAR  processing.  Amazingly,  it  is  as  simple,  if  not  simpler,  in  the  wavelet  case. 

Result  3.  Let  us  suppose  that  g  is  wavelet-domain  band-limited,  by  which  we  mean  that  there  exists  T>0 
such  that  for  all  s>T,  Wjg(u,v)  =  0  for  all  (u,v).  This  is  conceptually  similar  to  the  familiar  assumption 
of  band-limitedness.  Then  in  formula  (2)  above,  the  integration  with  respect  to  s  can  be  truncated  at  T 
instead  of  going  to  infinity,  and  it  is  possible  to  change  the  order  of  integration  and  integrate  with  respect 
to  5  first.  Integrating  with  respect  to  s  first,  but  making  the  change  of  variable  u=sx,  we  find  that: 

i^7?-LVe«**e)(f/(e))de,  (3) 

C/  0 

where: 

fio(v)  =  -y  uC  0(u)  du. 

That  is,  rather  than  needing  a  collection  of  returns,  it  suffices  to  filter  P^g  by  fie.  which  is  the  return  at 
angle  6  that  is  obtained  by  chirping  fie. 

It  can  be  shown  that  fie(v)  is  continuous  at  0,  and  with  proper  choice  of  wavelet  it  can  be  made  to 
decay  like  1  / v^.  Thus  it  is  entirely  reasonable  to  expect  that  fie  can  be  implemented  as  a  valid  chirp. 

We  see  that  once  the  convolution  has  been  performed,  the  only  other  computation  needed  is  integra¬ 
tion  with  respect  to  0.  Note  that  the  desired  chirp  fiafv)  depends  on  6,  so  that  potentially  different  signals 
are  needed  for  each  position  along  the  flight  path.  However,  if  the  wavelet  y  is  radially  symmetric,  then 
Ce  and  fie  will  be  independent  of  6.  Although  radially  symmetric  wavelets  arc  known  to  exist,  it  is  not 
known  if  these  wavelets  have  the  other  desirable  properties  required  to  make  fie  a  desirable  chirp  func¬ 
tion. 

Our  final  result  is  a  generalization,  and  returns  us  to  the  domain  of  Fourier  analysis,  as  opposed  to 
wavelet  analysis.  The  result  puts  the  previous  result  in  perspective,  showing  the  relationship  between 
existing  reconstruction  methods,  and  the  proposed  new  method. 

Result  4.  The  standard  reconstruction  formula  for  tomography,  and  implicitly  for  SAR  image  formation, 
is: 

=  ::^f’'(//'Feg)(xcose+>'sine)de. 

The  wavelet  form  of  the  same  equation  was  given  in  Result  1.  Here,  the  operator  H'  denotes  differentia¬ 
tion  followed  by  the  Hilbert  transfonn.  For  arbitrary  functions,  //'  cannot  be  implemented  as  a  convolu¬ 
tion  filter  in  the  spatial  domain,  and  must  be  achieved  by  multiplication  by  /p/  in  the  Fourier  domain. 
However,  if  we  suppose  that  /  is  band-limited  in  the  usual  sense,  then  application  of  H'  to  /  can  be 
achieved  by  convolution  with  any  function  q,  where  q(,x)  is  a  function  whose  Fourier  transform  equals  /p/ 
in  support  in  the  spectral  domain  of  /.  Tliat  is: 
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g(x,y)=~j^(P(.g*qKUmdQ,  (4) 

where  U  is  defined  as  before,  and  it  is  assumed  that  F^g  is  bandlimited  (as  a  function  of  «),  and  q  is  suit¬ 
ably  defined.  This  result  is  known  (for  example,  see  A.  Kak,  “Computerized  tomography  with  x-ray, 
emission,  an  ultrasound  sources,"  Proc.  IEEE,  September  1979),  but  not  widely  utilized.  We  sec,  how¬ 
ever,  that  Eqn.  (4)  has  the  same  form  as  Eqn.  (3),  and  implies  a  simple  method  for  SAR  image  formation, 
as  we  discuss  below.  However,  the  assumptions  are  different:  Eqn.  (3)  is  based  on  wavelct- 
bandlimitedness  of  g.  which  is  a  form  of  local-frequency  limitation,  whereas  Eqn.  (4)  assumes  the  usual 
bandlimitedness  of  Pog  for  all  6. 

We  now  discuss  the  implementation  of  foraiulas  of  the  form  (3)  and  (4)  for  SAR  image  formation. 
These  formulas  are  very  well-suited  to  SAR  for  the  following  reason.  All  existing  Radon  transform 
ipconstruction  technicpies  assume  that  the  values  of  the  Radon  transform  P^gix)  or  its  Fourier  transform 
Pegfp)  are  known  for  at  least  some  values  of  Q  and  x  or  p.  In  SAR,  however,  this  is  not  really  the  case; 
the  radar  signal  received  from  the  reflecting  object  is  not  P^g  but  PQg*S,  where  S  is  is  the  transmitted 
signal.  Using  a  filter,  it  is  possible  to  reconstruct  P^gip)  for  some  values  of  0  and  p  from  these  values, 
which  can  then  be  rerasterized  to  a  coordinate  grid,  in  order  to  obtain  g(^,Ti),  and  then  g  {x,y).  Clearly, 
considerable  processing  is  required  on  the  return  signals,  and  the  data  PQg*S  are  data  that  are  used  for  the 
reconstruction,  and  not  direcUy  related  to  the  image  data  in  an  untransformed  state.  With  our  methods, 
however,  convolution  of  Pgg  with  another  function,  which  we  may  call  the  convolving  function,  is  the 
first  step  in  the  reconstruction  technique.  We  therefore  recommend  implementation  of  any  of  these  three 
methods  by  using  its  convolving  function  as  the  transmitted  signal.  In  this  way,  the  convolution  of  P^g 


Figure  4.  An  outline  of  the  compulations  required  for  each  point  of  data  collection  in  the  SAR  image  for¬ 
mation  method  suggested  by  Eqns.  (3)  and  (4).  For  each  value  of  0.  the  return  signal  forms  a  line  of  data, 
which  is  spread  through  an  array  of  accumulators  in  a  direction  perpendicular  to  0. 
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first  step  in  the  reconstruction  technique.  We  therefore  recommend  implementation  of  any  of  these  three 
methods  by  using  its  convolving  function  as  the  transmitted  signal.  In  this  way.  the  convolution  of 
with  the  transmitted  signal  is  no  longer  a  nuisance  that  must  be  reversed,  but  the  first  step  in  the  solution! 

Accordingly,  the  method  for  reconstruction  is  as  follows.  At  each  point  along  the  llight  path  (see 
Figure  4),  a  convolution  function  h^{t)  is  sent,  and  the  return  signal  /e(u)  is  recorded.  For  any  given 
(jc,y),  the  contribution  of/e  toward  the  reconstructed  value  of  g(:ic,y)  depends  on /0(xcosGi>sin0).  That 
is,  we  have  the  formula: 

=  J  /e(Acos9+ysin0)d0. 

Note  that  the  contribution  for  any  given  /efu)  is  constant  along  a  line  perpendicular  to  the  direction 
(cos0,sin0).  Thus,  for  each  appropnate  value  of  «,  an  accumulator  along  a  line  perpendicular  to  the  direc¬ 
tion  0  is  incremented  by  an  amount  /efu).  This  can  be  done  in  parallel  for  different  values  of  u. 

To  summarize,  the  suggested  algorithm  is  as  follows: 

Initialize  g  (jc,y)=0  for  all  (x,y ). 

For  all  0  in  the  range  of  available  angles,  do: 

•  Obtain  a  return  signal  /efu)  using  a  chirp 

•  For  each  value  of  u,  do 

Use  Bresenham’s  algorithm  to  walk  along  a  line  in  a  direction  pcrpcrKlicular 
to  0,  (the  line  U  (8)  =  u,  and  increment  an  each  accumulator  by 
a  portion  of /efu). 

The  portion  should  be  computed  using  standard  anti-aliasing  methods. 

We  see  that  standard  computer  graphics  techniques  can  be  used  to  form  the  computations.  The  algo¬ 
rithms  are  easily  parallelize.  The  total  complexity  for  each  point  of  0  is  proportional  to  the  number  of 
pixels  in  the  image.  With  appropriate  hardware,  the  accumulation  could  take  place  during  the  collection 
of  data,  so  that  the  final  scene  would  be  available  by  the  time  the  final  data  point  is  collected. 
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